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Abstract 

In this paper, we study the macroscopic limit of a new model of collec- 
tive displacement. The model, called PTWA, is a combination of the Vicsek 
alignment model [H] and the Persistent Turning Walker (PTW) model of 
motion by curvature control |21|I24|. The PTW model was designed to fit 
measured trajectories of individual fish [23]. The PTWA model (Persistent 
Turning Walker with Alignment) describes the displacements of agents which 
modify their curvature in order to align with their neighbors. The derivation 
of its macroscopic limit uses the non-classical notion of generalized collisional 
invariant introduced in [20J. The macroscopic limit of the PTWA model 
involves two physical quantities, the density and the mean velocity of indi- 
viduals. It is a system of hyperbolic type but is non-conservative due to a 
geometric constraint on the velocity. This system has the same form as the 
macroscopic limit of the Vicsek model [20j (the 'Vicsek hydrodynamics') but 
for the expression of the model coefficients. The numerical computations 
show that the numerical values of the coefficients are very close. The 'Vicsek 
Hydrodynamic model' appears in this way as a more generic macroscopic 
model of swarming behavior as originally anticipated. 
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1 Introduction 

Modeling swarming behavior has attracted a lot of attention in the recent years. 
To model a flock of birds a school of fish [TBI [221 ESI H2] or the displacement 
of ants [mEnjHD], a key question is to understand how to relate the collective 
behavior of large groups of agents to simple individual mechanisms [TJUS]. From 
a mathematical point of view, this question takes the form of the derivation of 
macroscopic equations from individual based models [31 HI [TRESES]. This paper 
is devoted to the derivation of a macroscopic model for a new type of model of 
collective behavior where agents control their motion by changing the curvature 
of their trajectory. This model has been shown to provide the best fit of fish 
trajectories [21]. 

Among models of collective displacements, the so-called Vicsek model has re- 
ceived a particular attention [TS1HI]- This model describes the tendency of individ- 
uals to align with their congeners. Many features of this model have been studied 
such as the existence of a critical point [T21HT] . the long time behavior [T21I31] or 
the derivation of a continuum model [H12D]- Due to its simplicity, several exten- 
sions or modifications of this model have been proposed, such as the Cucker-Smale 
model [HI El CEO [23 [2B]- There is also a variety of models which add an attraction 
and a repulsion rule to the Vicsek model [T21I22]- However, the Vicsek model has 
been proposed on phenomenological bases. By contrast, the experiments of [21] 
have shown that the Persistent Turning Walker (PTW) model provides the best 
fit to individual fish trajectories. In the PTW model, the individual controls its 
motion by acting on the curvature of its trajectory instead of acting on its velocity. 
However, in its version of [2DI21], the PTW model only describes the evolution of 
a single individual. The model does not take into account the interactions between 
congeners. 

In the present work, interactions between individuals are introduced in the PTW 
model by means of an alignment rule, like in the Vicsek model. The resulting 
model, called PTWA (Persistent Turning Walker with Alignment) describe how 
each individual is influenced by the average velocity of its surrounding neighbors. In 
the framework of the PTW model where individuals control their motion by acting 
of the curvature of their trajectory, this influence must lead to a modification of 
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this curvature. This contrasts with the Vicsek model, where particles are directly 
modifying their velocity as a result of the interaction. 

The PTWA model is based on the assumption that the subjects use the time 
derivative of their trajectory curvature (or of their acceleration) as a control variable 
for planning their movement. Such models are not commonplace in the literature. 
Their first occurrence is, to the best or our knowledge, in [2T|I2^| . The present work 
is the first one in which interaction among the agents is taken into account within 
this kind of models (see also [25]). We note that [3B] introduces the acceleration 
of neighbors in the rule updating the subjects' velocities in a variant of the Vicsek 
Individual-Based model jH] but motion planning is eventually made by updating 
the velocity and not the acceleration. 

Once the PTWA model is set up, the main task of the present paper is to derive 
the macroscopic limit of this new model. This macroscopic limit is intended to 
provide a simplified description of the system at large scales. The major problem 
for this derivation is that there is nothing like momentum or energy conservation 
in the PTWA model. Such conservation laws are the corner stone of the classical 
theory of macroscopic limits in kinetic theory [IT1[TD]. Indeed, as a consequence of 
this absence of conservation, the dimension of the manifold of local equilibria in 
the PTWA model is larger than the dimension of the space of collisional invariants. 
Conservation laws are therefore missing for providing a closed set of equations for 
the macroscopic evolution of the parameters of the local equilibria. To overcome this 
problem, we use the notion of generalized collisional invariant introduced in [20] . 
Thanks to this new notion, a closed set of macroscopic equation for the PTWA 
model can be derived. 

The macroscopic model consists of a conservation equation for the local particle 
density and an evolution equation for the average velocity. The latter is constrained 
to be of unit norm. The resulting system is a non-conservative hyperbolic which 
shows similarities but also striking differences to the Euler system of gas dynam- 
ics. It has also the same form as the previously derived macroscopic limit of the 
Vicsek model (also referred to as the 'Vicsek Hydrodynamic model') in [2D], but 
for the expression of the model coefficients. At the end of the paper, we propose 
a numerical method to compute the generalized collisional invariant out of which 
the coefficients of the macroscopic model are derived. The similarity between the 
'Vicsek hydrodynamics' and the 'PTWA hydrodynamics' can be better understood 
by considering the relations between the microscopic models. Indeed, the Vicsek 
model can be seen as a special limit of the PTWA model in a well-suited asymptotic 
limit. Work is in progress to establish this connexion firmly. 

The inclusion of the alignment rule in the PTW model changes drastically the 
large scale dynamics of the system. Without this alignment rule, the PTW model 
exhibits a diffusive behavior at large scales [TD1[2T]. By contrast, when the align- 
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ment rule is included, the model becomes of hyperbolic type. Therefore, the local 
alignment rule added to the PTW model generates convection at the macroscopic 
scale. 

Since the addition of the alignment rule modifies drastically the dynamics of the 
PTW model, it is also interesting to study the large scale effects of other types of 
local rules such as attraction-repulsion. The goal is to find a common framework 
for the large scale dynamics of a large class of swarming models. Currently, there 
exist a profusion of individual based models, especially for fish behavior (see for 
a short review). In a macroscopic model, only the gross features of the microscopic 
model remain. Therefore, the derivation of macroscopic models may be a tool to 
better capture the common features and differences between these different types 
of swarming models. 

The outline of the paper is as follows: in section 2, we introduce the PTWA 
model and the main result is stated. Section 3 is devoted to the proof of the 
derivation of the macroscopic limit of the PTWA model. In section 4, we study 
some properties of the so-obtained macroscopic model and we numerically estimate 
the involved coefficients. Finally, in section 5, we draw a conclusion of this work. 

2 Presentation of the model and main result 
2.1 The individual based model 

The starting point is a model in which alignment interaction between agents is in- 
troduced inside the Persistent Turning Walker model (PTW) [2TJI23]- The PTW 
model is a model for individual displacements which has been derived to fit exper- 
imentally observed trajectories of fish. It supposes that individuals control their 
motion by acting of the curvature of their trajectory. To make it a realistic model 
for collective displacements, the PTW model must be enriched by introducing inter- 
individual interactions. Indeed, one of the main features of collective motion such 
as those observed in animal populations (fish schools, mammalian herds, etc.) is 
the ability of individuals to coordinate with each other. Observations suggest that 
trend to alignment is an important component of this interaction and leads to a 
powerful coordination-building mechanism by synchronizing the agent's velocities 
one to each other. One of the simplest models of alignment interaction is the Vicsek 
model [JT]. This time-discrete model supposes that individuals move at constant 
speed and align to the average velocity of their neighbors (up to a certain stochas- 
tic uncertainty) at each time step. A time-continuous version of this dynamics has 
been derived in [20J. 

In order to combine the PTW displacement model and the Vicsek alignment 
interaction model (in the time-continuous framework of [2D]), we propose the fol- 
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lowing model further referred to ad the PTWA model (PTW model with alignment): 
among a population of N agents, the motion of the i th individual is given by 

£ - <*«.), 

ddi 
dt 

dtii = a(vKi — Ki) dt + b dB\, 

with 

^ = r(9i) x Hi (2.4) 

and 

a = 4r , J* = c£ r(0j), (2.5) 

where x = (xi,x 2 ) G M 2 is the position of the individual, r(^) = (cos# i; sin^) 
is the direction of its velocity vector, with the angle 6\ G (— 7r, 7r] measured from 
the Xi direction, k ; e 1 is the curvature of its trajectory and B\ is a standard 
Brownian motion (with B\ independent of B\ for i ^ j). The magnitude of the 
velocity is constant and denoted by c > 0. The constant a is a relaxation frequency 
and b quantifies the intensity of the random perturbation of the curvature. The 
vector Qi is the mean direction of the neighbors of the i th individual (defined as 
the individuals j which are at a distance less than R from x$, R > being the 
perception distance of the individuals, supposed given). 

The trend to alignment is modeled by the relaxation term of (12.31) (in factor of a). 
It describes the relaxation of the trajectory curvature to the target curvature Ki. Ki 
is computed by taking the cross producl0of the direction of the individual f{Bi) and 
the mean direction of its neighbors VLi. vKi is the trajectory curvature the individual 
must achieve in order to align to its neighbors. It increases with increasing difference 
between the individual's velocity and the target velocity, v is the typical value of 
the individuals' trajectory curvature and can be seen as the 'comfort' curvature. 
The larger v is, the faster alignment occurs. The second term of (12.31) (in factor of 
b) is a random term which describes the tendency of individuals to desynchronize to 
their neighbors in order for instance, to explore their environment. At equilibrium, 
these two antagonist effects lead to a stationary distribution of curvatures which is 
the building block of the construction of the macroscopic model. 

We illustrate this model in figure [U In the left figure, a fish is represented 
turning to the left. However, its neighbors are moving towards the other direction 
(Q is pointing to the right). Then the fish is going to adjust its curvature in 

^or two-dimensional vectors a = (01,02), b — (61,62), the cross product a x b is the scalar 
0162 — 0261. 
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(2.2) 
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order to move towards the same direction as Q (right figure). The adjustment of its 
curvature requires a certain time of order 1/a, after which the curvature k is close to 
vK. Therefore, in this model, there is a time delay between the current acceleration 
of the fish (k) and its desired acceleration (vk). In most models describing animal 
behavior, the dynamics is inspired by Newton's second law: the acceleration of an 
individual is equal to a force term which incorporates all information about the 
environment. In the present model, individuals need a certain time to adjust their 
acceleration. This rule can be seen as a modification of Newton's second law saying 
that the force is proportional to the time derivative of the acceleration rather than 
to the acceleration itself. 



Figure 1: Illustration of the model ( 12. II) — (12.30 . On the left figure, a fish is turning 
to the left, while its neighbors are moving to the right (Q). After a certain time of 
order 1/a, the fish adjusts its curvature in order to align its velocity with Q (right 
figure). 

Our goal is the study of model (12. 1 1) — (1 2 . 3 j) at large time and space scales. For 
this purpose, it is convenient to introduce scaled variables. We use xq = v~ x as 
space unit, to = (ci>) _1 as time unit, k = Xq 1 = v as curvature unit. We introduce 
the dimensionless time, space and curvature as t' = t/to, x' = x/x and k' = k/ko 
and for simplicity we omit the primes in the discussion below. In scaled variables, 
the PTWA model is given by (for the i th individual) : 



r 





(2.6) 



dt 



(2.7) 
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with Ki defined by equation (l2.4IUI2.5p (c being replaced by 1) and A, a given by: 

a 2 b 2 

A = — , a - 



cv 2cv 3 
2.2 Main result 

A first step consists in providing a mean-field description of the PTWA dynam- 
ics. Introducing the probability density function of fish /(t, x, 0, k), we will prove 
formally that the PTWA model (I2.6p - (I2.8I) leads to the following equation for /: 

d t f + m • V x / + Kdgf + \d K [(« -«)/]= a 2 d 2 J, (2.9) 

with 

K = f(fl)xfl(x) (2.10) 

and 

'|x-y|<R,0,/e 



fi ( x ) = T7rTi > J ( x )=/ f(9)f(y,e,K)dyded K . (2.11) 



The main concern of this paper is the study of the so-called hydrodynamic limit 
of the mean-field model (I2.9p . With this aim, we perform a new rescaling and 
introduce the macroscopic variables t and x: 

i=et , x = ex, (2.12) 

with e > a small number representing the ratio between the microscopic and the 
macroscopic time and space scales. In this paper, we give a formal proof that the 
density distribution of individuals in these new variables f e (t,5t,9, k) converges in 
the limit e — > to the solutions of a hydrodynamic like model. More precisely, the 
theorem reads (dropping the tildes for simplicity): 

Theorem 2.1 In the limit e — > 0, the distribution f e converges to an equilibrium: 

f P Mu(e)M(K) 

with and Af (resp.) a Von Mises distribution and a Gaussian distribution 
defined at Ii3.12\) and \3. 1(J\) . The density p = p(x,t) and the direction of the flux 
Q = Q(x, t) satisfy the following system: 

d t p + cxVx ■ (fjQ) =0, 

p(d t vt + c 2 (n ■ v x )ft) + f^(id - n <g> n)v xP = o, (2 ' 13) 

where C\ and c 2 are two positive constants defined later on at 113.29]) k3~c 
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The so-obtained macroscopic model (I2.13P has the same form as that derived 
from the Vicsek model [2U|. Indeed, the two models only differ by the values of 
their coefficients. This model is a hyperbolic system which bears some similarities 
with the Euler system of isothermal compressible gases. There are however some 
striking differences. First, the convection speed of the density p is different from the 
convection speed of the velocity Q (c\ ^ in general). Moreover, the velocity Q is 
a unit vector and therefore it satisfies the constraint |Q| = 1. This explains why the 
pressure term is premultiplied by the matrix (Id — Q ® Q). This projection matrix 
guarantees that the resulting vector is orthogonal to Q. Consequently, the constraint 

= 1 is preserved by the dynamics. However, the projection matrix leads to a 
non- conservative model which cannot be put in conservative form. This intrinsic 
non-conservation feature is the macroscopic counterpart of the lack of momentum 
conservation at the microscopic level (see below). 

The modification of the PTW model leading to the PTWA model has drastically 
changed the nature of the macroscopic model. Indeed, the macroscopic limit of the 
PTW model without the incorporation of the interactions is of diffusive nature 
PHIEI]- By contrast, that of the PTWA model is of hyperbolic type. Indeed, the 
scaling (12.121) is of hydrodynamic type, the macroscopic time and space scales being 
of the same order of magnitude. By contrast, a diffusive scaling would have required 
t = e 2 t instead (see [T01I2T]). 

The similarity with the 'Vicsek Hydrodynamics' also confirms that that the 
chosen interaction rule generates alignment since the PTWA model has the same 
macroscopic limit as the Vicsek model. At the microscopic scale, the PTWA and 
Vicsek models look rather different, whereas, at the macroscopic scale, they are 
similar. This is an example of how the derivation of macroscopic model can be used 
as a tool to reduce and unify different types of swarming models in classes leading 
to similar macroscopic models. 

3 Derivation of a macroscopic model 
3.1 Mean field equation 

In this section, we briefly summarize the first step of the derivation of the macro- 
scopic model, namely the derivation of the intermediate mean-field equation (12. 9p 
from the particle dynamics (I2.6l) - (l2.8p . In order to derive this mean field equation, 
we start by looking at the system without the white noise dB\ for a large (but fixed) 
number of individuals N. In this case, the system reduces to a coupled system of 
ordinary differential equations. We denote by {X i {t),@i{t),K i (t)']i = i,,,N the solu- 
tion of this system on a given time interval. Following the standard methodology 
(see e.g. the text book [57]). we introduce the so-called empirical distribution f N 
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given by: 



N 



N ^ 



Xiit) ® <W) 



5 



Ki(t)- 



(3.1) 



i=l 



We can easily check that this density distribution satisfies the following equation 
(weakly) : 



d t F +r(0)-V x r +Kd e f +^d K 



A' 



0. 



with 
and 



N , 



K N = t{6) x Q (x) 



3 N ( 



E 



I V /I j, |x-x 3 -|<J? 

The term J N can be expressed using the empirical distribution / : 

r(Q)f N (y,0, «) dyd0d/c. 



J^x) 



AT 



'|x-y|<_R, 0,k 

Then it is clear that the formal limit iV — > oo of / satisfies the following equation: 



d t f + f{9) ■ V x / + Kd e f + A<9 K (7c - n)f 



0. 



with 7c given by (j2TTUl) . (j2TTT|) . 

When the white noise is added, the situation is more complicated. At the particle 
level fl2.6p - fl2.8l) . the system becomes a coupled system of stochastic differential 
equations. This implies that the empirical distribution / given by (13.11) becomes 
a stochastic measure. In this case, formal considerations suggest that, in the limit 
A" — > oo, the distribution function / satisfies the following Fokker-Planck equation 
(|2.9p with 7c given by fl2.10l) . fl2.11l) . For related questions, we refer the reader to 

MM]. 
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In order to derive a macroscopic equation from the mean- field equation fl2.9p - fl2.lll) . 
we use the hydrodynamic scaling. With this aim, we introduce the macroscopic 
variables t and x defined by (12. 12ft . In the rescaled variables, the distribution 
function (denoted by f £ ) is given by f £ (t, x, 9, /c) = ^ f(t, x, 9, k). After omitting 
the tildes, it satisfies the following equation: 



e(d t f + r(9) ■ V x / £ ) + Kd e f £ + Xd K (Tc^ -n)f = a 2 d 2 J 



(3.2) 



with 



K 



f{9) x n £ 
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and 

n£ ( x ) = ni7^T > J£ ( x )=/ r(e)f(y^^)dydedK. (3.3) 

|J V X JI J\x-y\<eR,8,n 

We note that the expression (13.31) of J e supposes that the radius of interaction be- 
tween the individuals is tied to the microscopic scale. This assumption translates 
the fact that in most biological system, each individual has only access to infor- 
mation about its close neighborhood. Thanks to this assumption, we can replace 
the expression of Q by a local expression. This is precisely stated in the following 
lemma, the proof of which is obvious and omitted. 



Lemma 3.1 We have the expansion: 



where 



ft" = Qfe + 0(V 



Q f ^) = A^L and j e (x)=/ T{6)f £ {^0^)d6dK. 
|J l x JI j o,k 

Finally, we can simplify f)3.2p using the equality: 

f (0) x n = sm(6 - 6) 

with 6 such that: 

With these notations, equation ( 13. 2h can be written as: 

e(d t f + f{9) ■ V x / £ ) = Q(f) + 0(e 2 ) (3.4) 



with the operator Q (below referred to as the 'collision operator') defined by: 

Q(f) = -Kd f - X sm(6 - 6) dj + \d K ( K f) + a 2 d 2 J, (3.5) 
where f{6) = fi/(x) defined as: 

n f (*) = 0L , j(x)=/ T(6)f(x,6,K)d9dK. (3.6) 

In the sequel, we will drop the 0(e 2 ) remainder which has no influence in the final 
result. 
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3.3 Study of the collision operator 
3.3.1 Equilibria 

In order to study the limit e — > of the solution f £ of (13. 4ft . we first have to 
determine the equilibria of the operator Q defined by ( 13.51) . With this aim, we 
notice that Q can be decomposed as a sum of a formally skew-adjoint operator 
and of a formally self-adjoint operator. For the skew-adjoint part, we introduce the 
function: 

K 2 

H(6,k) = -\cos9 + — 

and we adopt the convention that for any function h(9, k): 

h n (9,K) = hg(9,K) = h(9-9,K), (3.7) 

with t(9) = Q. Using these notations, for any smooth function /, the skew-adjoint 
part of Q can be written as: 

- Kdgf - A sin(£ - 9) dj = d d H w dj - d K H- e d e f = (3.8) 

using the Poisson Bracket formalism {-, -}(0 )/t ) in the (9, k) space. Therefore, any 
function of the form g(Hj) satisfies {Hg, g(Hg)} = 0. On the other hand, the 
self-adjoint part of Q satisfies: 



Xd K (Kf) + a 2 d 2 J = a 2 d K (Afd K Uf\\ 



(3.9) 



with M the Gaussian distribution with zero mean and variance a 2 /X: 



^ - ^ (-^) • < 3 - 10 > 

In particular, the Gaussian Af is in the kernel of the self-adjoint part of Q. We 
combine our two previous observations to define the function: 

/i(#,K) = Cexp^-A#) =Cexp(-^ (^-\cos9j\ , (3.11) 

where C is the normalization constant such that j, e R ^ fi(9, k) d9 dn = 1. This nor- 
malization constant is explicitly given below. The translates fig of /i in the sense 
of definition (13.71) are of the form g(Hg) and are Gaussian distributions in k with 
variance a / VX It follows from a simple computation that /ig is an equilibrium for 
Q (i.e. Qifig) = 0), for all real values of 9. 
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To simplify the analysis, we introduce the Von Mises distribution M.: 

M{6) = C exp (^cose\ , (3.12) 

where Co = (27rJo(^))~ 1 is the normalization constant such that J B Ai(9)d9 = 1 
(with I the modified Bessel function of order 0). Therefore, p can be written as 
the product of M. given by (I3.12|) and M given by (I3.1U|) : 

h{B,k) =M{9)M{k), (3.13) 



and the normalization constant C is given by C = C y A/(27ra 2 ). We summarize 
our analysis of Q in the following proposition. 



Proposition 3.2 i) The operator Q satisfies: 



[ Q(f) ^- dOdn = -a 2 f 

J9,K Lin J9,K 



d K 



'r 

.AC 



2 



dddn < 0, (3.14) 



p B Je,K 

with p defined by Ii3.13i) and 9 such that f{6) = Qf with Qf defined in ( 13.61) . 

ii) The equilibria of Q (i.e. the functions f(9,n) > such that Q{f) = 0) form 
a two-dimensional manifold £ given by: 

E = {pp- e | P eR + ,0e (-tt, tt]}, (3.15) 

where p is the total mass and 9 the direction of the flux of p p^. 
Proof, (i) Combining (13.81) and f)3.9p . we find: 

Q(f) = {H- 6 , /} + a 2 d K Ufd K (LYj . (3.16) 

Using (13.111) . the fact that the Poisson bracket with / is a derivation and is a 
skew-adjoint operator, we find: 



/ {HnJ}-d9dK = [ ^e-^ H e{e^ H e,f}td6dK 
Je,K p Je,n A p 

= T-h f {^ H °,f}fd&dK = 0. 
a o Je,n 



Then, using the formulation of Q in ( 13.161) . we easily deduce the equality (I3.14p by 
applying Green's formula. 
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^ii) If / is an equilibrium for Q (i.e. Q(f) = 0) using the equality f)3.14p we have: 



2 



dOdn = 0, 



which means that / is proportional to M as a function of n. Therefore, we can 
write: 

f(8,K)=cp(8)Af(K). 

Using again that / is an equilibrium, we have: 

-Ktp'(9) + A sm(9 - B)^ip(B) = 0, for all k. 

cr 

Solving this differential equation leads to if (9) = CMq{9) with M. given by (13.121) . 
This yields / = KAigAf with K > a constant which proves that / is of the form 
/ = p/i0 o , with p > and 9 £ (— n, ir]. 

Reciprocally, we show that a function of the form / = pfig with p > and 9$ £ 
(— 7r, it] is an equilibrium. For this purpose, the only thing to show is that the 
associated Qf = r{9) is such that 9 = 9q. We compute 



j/ = / ppe r(9)d9dK 

= p k«) M{K)C(i exp i$ co < d - 9 ^) ( Zl ) d6dK - 

Then, by the change of variables <fr = 9 — 9q and using oddness considerations, we 
obtain 



if = pCo exp ( 



2 COS</> COSpCkp T{V ) = p " 2 T{V ), 



where l\ is the modified Bessel function of order 1. Remembering that Qf — j//|j/|, 
we deduce that 9 = ±9 , with the sign being that of h(^)/Io(^)- A simple 
inspection of the integral giving Ii shows that this sign is positive and that 9 = 9 , 
which ends the proof. □ 



3.3.2 Generalized collisional invariant 

The next step to determine the hydrodynamic limit of f £ ( 13.41) is to look at the 
collision invariants of the operator Q, i.e. the functions ip which satisfy: 

/ Q(f)ifjd9dK = 0, for all/. 

J0,K 
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Clearly, ij; = 1 is collisional invariant. But there is no other obvious collisional 
invariant. However, since the equilibria of Q ( 13.151) form a two dimensional space, 
we need two conserved quantities to derive a macroscopic model. To overcome this 
problem, we use the notion of generalized collisional invariant developed in [2D]. 

In this paper, we use a slightly different definition from [20J. Indeed, the result 
of [2D] was slightly incorrect and the present definition is designed to make the 
statement correct. We first introduce the following definition: 

Definition 1 For a given Q £ S 1 and a given distribution function f(8,K), we 
define the 'extended' collision operator Qn(f) by: 



Gn(/) = {Hn, /} + a 2 d K (W (Ljj 



where we recall the notation 113. 7\) . 
Obviously, we have 

Q(f) = Qa f (f), (3.17) 

recalling the definition (13. 6p of Qf. For fixed Q, the operator Qn(f) is linear. We 
now define a Generalized Collision Invariant. 

Definition 2 For a given unit vector Q £ S 1 , a function ipn is called a Generalized 
Collisional Invariant (GCI) if it satisfies: 

I Qn(f)?Pnd9dK = 0, for all f such that Q f = ±Q, (3.18) 

Using definition (13. IS)) with Qf = Q and (j3.17]) : we note that if ipn is a GCI, it 
satisfies 

Q(f) ^ nf d6dK = 0. 



ie,K 

This property is crucial for the establishment of the hydrodynamic limit. 
For a given Q £ S 1 , the adjoint operator to Qn is given by: 

Q* u (i>) = ndei) + A sin(^ - 9) d K ^ - Xnd^ + a 2 c^, 

with 9 such that Q = r{9). This operator Qq enables us to find an explicit equation 
for the GCI ipQ as stated in the following lemma. 

Lemma 3.3 For a given unit vector Q £ S 1 , a function ipn is a generalized colli- 
sional invariant if and only if it there exists a constant (3 £ M such that: 

QMn) = Pm xa (3.19) 
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Proof. Let f(9, k) be such that Qf = ±fi. This is equivalent to saying that there 
exists a constant C 6 I such that j/ = CQ (see ( 13.61) for the definition of jf), or 
in other words, that jf x Q = 0. Now, if ip satisfies ( 13.191) . we have, for such a 
function /: 

/ Qn(f)^pd6dK = I f Q*(V0 d9dK 

JQ.n J9,k 

= (3 j ff{9) x ft dOdn = j3j / xfi = 0, 
and ^ is a GCI associated to Q. 

Reciprocally, if i/jq is a GCI associated to Q, we have: 

/ Qa(f)ipad6dK = 0= [ fQ* n (^n)d9dK 

for all f(9, k) such that j/ x Q = 0. We deduce that, for all /, 

j /X fi = =^ / fQ* {l {i/> sl )d6d K = 0. (3.20) 

The two expressions appearing in ( I3.20p are linear forms acting on /. By an el- 
ementary lemma [B], the one appearing in the right-hand side is proportional to 
the one appearing in the left-hand side, with a proportionality coefficient (3 G M.. 
Expressing this proportionality gives: 

f{Qh(^)-r{9)xQ)d9dK = 0, (3.21) 

for all / without any restriction. ( 13.211) yields ( 13. 19ft . which concludes the proof. □ 

It remains to prove the existence of GCI's, or, in other words, to prove the 
existence of solutions to equation (13.191) . With this aim, we use the Hilbert space 
L 2 equipped with the scalar product < ., . > fl defined by: 

Ll = {f(6,K)/ ( \f\ 2 fxd9dn <+oo}, 
<f,9>n= I fgi^dOdK. (3.22) 
Below, we will also use the notation: 

<g>^= f 9 (9,k) n(9,K)d6dK. (3.23) 

J9,K 

We define the hyperplane E: 

E = {f eLl(9,K) j I ff,d9dK = 0} 
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and the linear operator C: 



Cip = ndgijj 



A sin 6d K ifj — \K,d K ip + a 2 d 2 ip 



(3.24) 



with domain D(C) given by: 



D{C) 



{feLl/Cfe Lj}. 



We have the following lemma: 



Lemma 3.4 (i) Let \ G L 2 . A necessary condition for the existence of a solution 
ip G D(C) of problem 




= 0. 



(ii) For all \ G E, the problem A3.25\) has a unique solution ip in E. Then, all 
solutions to problem A3.25\) are of the form ip + K, with an arbitrary K G M. 

In Appendices Al and A2, we give two different proofs of the fact that (13.251) is 
uniquely solvable in E. The proof in appendix Al uses tools from functional analysis 
(see also [21J). The proof in appendix A2 uses probabilistic tools to analyze the 
stochastic equation associated to (I3.25P (see also [ID])- Here we only prove (i) and 
the last statement of (ii). 

Proof, (i) The formal adjoint C* of C is given by the expression (13. 5p of Q in which 
9 = 0. Therefore, from section 13.3.11 we have that £*(/i) = 0. Integrating (13. 25ft 
against /i and using Green's formula leads to the necessary condition J R x^ d9dn = 
0, i.e. to the fact that x must belong to E. 

The second part of (ii) amounts to showing that the null space of £ reduces to the 
constant functions. Indeed, it is straightforward to see that £(1) = 0. To prove 
that the constant functions are the only elements of the null space of C, we suppose 
that if) G D(C) such that Cip = 0. Using that < Cip, ip > M = 0, we find, using 
Green's formula: 



Therefore, i/j is independent of k. So we can write: ip(8,K) = Using again 

that £$ = 0, we find that $ is a constant. 




We refer to appendices Al or A2 for the existence part of point (ii). 



□ 



The following proposition completely determines the set of GCI's associated to a 
vector Q. 
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Proposition 3.5 For a given f2 £ § , tie set of the GCI's associated to f2 is a 
two dimensional vector space Cq = Span{l,^n} where ^ is given by: 

ifa(6,K) = if>(6 (3.26) 
with 9 such that r(9) = Q and ip is the unique solution of: 

Cip = -sm6, (3.27) 
belonging to the hyperplane E. Moreover, the function ip satisfies the property: 

^(-0,-«) = -<0(0,«). (3.28) 

Proof. We first note that (13.191) is a linear problem and that it is enough to solve 
it for (3 = 1. Simple calculations show that ipn is a solution to (I3.19P if and only if 
there exists a function ip such that i^q(9) = ip{9 — 9) with ip a solution of (13.271) . 
This shows (ET2SI) . 

To show the existence and uniqueness of a solution ip to (I3.27P in E, it is enough 
to check that the right-hand side of (13.271) belongs to E i.e. satisfies the compat- 
ibility condition J 9 n xiid9dK = 0. But this follows readily by oddness considera- 
tions. Moreover, noting that the operator L is invariant under the transformation 
(9, k) —> (—9, —k), (I3.28P follows from the uniqueness of the solution. 

Again, by the uniqueness in E and by the second part of Lemma 13.41 (ii), all 
solutions to (13.271) consist of linear combinations of ip and of a constant function. 
It follows that the set of GCI's associated to Q is the two-dimensional vector space 
spanned Cq = Span{l,ij)n}. This ends the proof. □ 

3.4 Limit e -> 

Since we know the equilibria and GCI's of the operator Q, we can give a formal 
proof of theorem 12.11 

Proof of Theorem 12.11 If we suppose that f e converges (weakly) to f° as e — > 
we first have: 

W°) = o, 

which means that f° is an equilibrium. Thanks to section 13.3. 1\ f° can be written 
as: 

f = fPMa>iW{K), 

with M. and J\f defined in ( 13.121) ( I3.10p . The mass p°(t, x) and the direction of the 
flux Q°(t, x) are the two remaining unknowns. 
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In order to find the system of equations which determines the evolution of p° 
and Q°, we first integrate (13 .4p with respect to (9, k). We find the mass conservation 
equation: 

<9*p £ + V x -j £ = 0, 

with 

■f = ( f(9)f £ d9dn. 

In the limit e — > 0, this gives: 
with the constant C\ given by: 



r I (—) 

ci = / cos9M(9)d9 = (3.29) 
J® lair?) 



Therefore we deduce that p° and Q° obey the following mass conservation equation: 

<V + c x V x • (p°fi°) = 0. 

In order to fully determine the evolution of p° and Q°, we need to find a second 
equation. For this purpose, we integrate (13. 4p against the generalized collisional 
invariant ipQe (I3.26p . with Q £ = Qfe. This leads to: 

/ (d t f + t(6) ■ V x f)^ dOdn = 0. 

In the limit e — > 0, we find : 

/ d t (p°M n oAf) ip u0 d6dn + / f(9) ■ V x (p°M n oN) tp n o dOdn = 0. (3.30) 
Je,K Je,K 

For clarity, we drop the exponent '0' and write (p, Q) for (p°, Q°) in the discussion 
below. Using polar coordinates for Q = f{9) = (cos 9, sin 6), elementary computa- 
tions show that: 

d t (pMo) + f(9)-V x (pMg) = d t pMe + pMe^sm(9-9)d t 9 

+ f{9) ■ (V xP Me + pM- e ^ sin(# - 9) vj) . 
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Therefore, equation (I3.30p leads to: 

d t p MjAf ipg dOdn 



+ J t{6) ■ (y^pMgMi^) dOdn 

+ — / t{6) ■ (p MqM sin(fl - 9) VJ ipg) dOdn = 0. 



This equation can be simplified using the symmetry satisfied by ip (I3.28p . We treat 
each term separately. First, we have: 

X 1 = I d t pM^ipgdedK 

= d t p [ M(6-lftAr(K)M6-6,K)d6dK = 0, (3.31) 
Je,K 

because M.{6)N(k) is an even function of the pair (9, k) and tp(9,n) is odd. For 
the second term, we use the change of unknowns 9' = 9 — 9 and get: 

X 2 = -^pd t 9 M(9')N(K)sm9'm',K)d9'dK 
cr Jo',k 

%pd t 9 lx , (3.32) 



a 2 



with 

7i=<sin^^> M (3.33) 
using the notation (I3.23p . For the third term, we find: 

= V x p- / T(9)MsAftkgd9dK 

J 8,K 

= V x p- / T(9 + 9)M(9W(K)ip(9,K)d9dK 

— f ( cos 9 cos 9 — sin 9 sin 9 \ . .... . ,, N , . ,„ , 

= W^sm^ + cos^ )Mm{«)m«)MdK. 

Once again, using the symmetry satisfied by ip, we find: 

„ / -sin# \ 
X 3 = 7iV xP - cosfl ' 
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with 71 defined in (I3.33p . If we denote by r(6 l )- L = Q 1 - the orthogonal vector to 
f(5): 

y cos0 J 

we finally get: 

X 3 =7iV x p-r(^) ± - (3-34) 

For the last term, we have: 

\2 , 

X 4 = — p V x • / T(8)Mg(8)J\f(K) sin(0 - 0)^(0, «) d^d/c 

\2 . 

= -rpV x 0- / t(6 + 6)Mo(0W(k) sm6 ip(8,K)d9dK 
cr Je,K 

= ^7 2 pV x 0-f(0), (3.35) 
cr 

with 

72 =< cos 6* sin0?/; > M . 
Combining fEPTD . fl3T32D . fECT) and (I3T35D yields: 

Ti^P^ + TiVx/o-W + 72^/° V x^-^) = 0. (3.36) 

cr cr 

Using again the unit vector Q = r(0), elementary computations show that: 

d t n = d t 9n ± and (tt ■ V x )fi = (fi x ® 
Therefore, multiplying equation (I3.36P by Q 1 - leads to: 

p5 t Q + — (Vxp-n 1 )^ 1 + — pin- V x )fi = 0. 
A 2 7i 

This finally leads to: 

a 2 

p d t n + c 2 p(n- v x )fi + — (id - n ® n) v x p = o, (3.37) 



with 

72 < sin cos ip >„ 
c 2 = — = ~ : 7TT~r » (3-38) 
7! < sin V >m 

which end the proof. □ 
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4 Properties of the macroscopic system 
4.1 Hyperbolicity 

The macroscopic system (12.131) arising from the PTWA dynamics has the same form 
as the system found in [20] for the macroscopic limit of the Vicsek model. Indeed, 
if we define the diffusion coefficient d as: 

a 2 

2 

then the coefficient c\ given by ( 13. 29ft and the coefficient % in front of the pressure 
term in (13.371) are exactly the same in the two systems. Only the coefficient c 2 given 
by (13.381) differs from that of [20]. Thus, the study of the hyperbolicity of system 
(I2.13P is completely similar to the one conducted for the Vicsek model in [2"0"II3"3"] . 
We briefly summarize the analysis here. Using the geometric constraint \Q\ = 1, we 
can parametrize the direction of the flux Q in polar coordinates: Q = (cos 9, sin 9) 
with 9 g] — 7r, 7r]. In order to look at the wave propagating in the ^-direction, we 
suppose that p and Q are independent of y. Therefore, under this assumption, the 
system ( I2.13P reduces to: 



dtp + C\d x (p cos 9) = 0, 
a 2 sin 9 
A 2 " 

The characteristic velocities of this system are given by: 



d t 9 + c 2 cos 9d x 9 - -r^—dxP 



1 

7= 2 



/ Q 2 

C\ + c 2 ) cos 6 ± \ (ci — c 2 ) 2 cos 2 9 + 4c! — sin 2 9 



A- 



The system is therefore hyperbolic since the characteristic velocities are real. 



4.2 Numerical computations of ip 

In order to compute the macroscopic coefficient c 2 (I3.38p . we first need to calculate 
the generalized collisional invariant ip (I3.26p . With this aim, we introduce a weak 
formulation of the equation satisfied by ip. In the Hilbert space -^(S 1 x M), the 
function ip satisfies: 

< Ctjj, if > M = - < sin^, > M , V^GLj, (4.1) 

where the scalar product < ., . > M is defined in (I3.22p and the operator £ in (13.241) . 
To approximate the solution ip numerically, we use a Galerkin method. It consists 
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in solving the weak formulation (14. ip for all the functions (f in a subspace V of 
of finite dimension. To construct such a subspace V, we use a Hilbert basis of L^. 
For this purpose, we consider the following functions: 



2tt M(6) 



PJk) 



K 



where Ai is defined in ( I3.12p and H n is the n th Hermite polynomial. We can easily 
prove that the family {(p m P n }m,n>o is a Hilbert basis of IA. Then, for any odd 
positive integers m and any positive integer n, we define the vector space V m ,n- 

V m , n = Span{(/?jPfe / \j\ < m , < k <n}. 

The Galerkin method consists in finding ^ m>n G V m>n such that equation fl4.ll) is 
satisfied for every <p G V m , n '- 



ra.ni 

We can decompose ip as: 



< sin 6, ip > M , V<£> G Kn, re . 



\j\<m,0<k<n 



(4.2) 
(4.3) 



where Cj 1 are complex coefficients given by: 

Cj =< i) m ,n, <Pj Pk >fi ■ 

We store the coefficients {C|}i,i< m .n<Kn in a matrix X such that: 



y i J lil<m,0<fc<n 

X{j,k) 



(4.4) 



We call the matrix X the matrix representation of ip m ,n in Kn,n- We want to 
transform the problem satisfied by ij) mn (14. ip into a matrix equation for X. With 
this aim, we define several matrices. 



Definition 3 We define the matrices L_i and L +1 by: 

1 



L-i 





1 



and the diagonal matrices: 

D 2 -- 



J +i 



diag(-m, ...,-1,0, 1, 
diag(0, 1, 2, ... n). 



m) 



1 





(4.5) 
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Using the matrices defined above, we can convert the equation satisfied by tp Ti 
into a matrix equation for X. 



Proposition 4.1 Let ij} mn G V m) n the solution of I j4.1)) in V m>n . Its matrix repre- 
sentation X = {Cj }\j\<m,o<k<n m the Hilbert basis {<p m P n } satishes: 



with 



P-lM-lXN-l + f3 2 M 2 XN 2 - \XD 2 = B (4.6) 



M x = D 1 , JVi = y/L\L^ + L +iy /D;, ( 4 - 7 ) 



M 2 = L_i-L +1 , N 2 = y/I^L-! - L +iy /D 2 
and B the matrix representation of — sin 9 in V mtTl given by: 



where L is the modified Bessel function of order j . 



otherwise, 



Since the demonstration of proposition 14.11 is only a matter of computations, we 
postpone the proof to appendix B. To solve (14. 6p . we transform the linear equation 
( 14. 6 p into a linear system that we invert numerically. This eventually allows us to 
construct ip m ,n using (14.31) . 

On figure [2] (left), we display an example of an approximate solution il) m ,n of the 
GCI ij) for A = 1 and a — 1. We also estimate Cip m ,n numerically using a finite 
difference method (figure [2J right). The figure clearly suggests that Cip m ^ n is close 
to — sin 9, providing a qualitative check of the accuracy of the computation. To 
make this assessment more quantitative, we compute the residual \L^) m<n + sin^loo 
for different values of (A, a) on figure [31 As we can see, the residual gets larger 
when a increases and gets smaller when A increases. 



4.3 Computation of the coefficient C2 

Once we have computed the generalized collisional invariant ip, we can calculate 
the coefficient c 2 using (13.381) . On figure @J we fix the the parameter A = 1 and we 
compute the value of c 2 for different values of a (we still use m = 30 and n = 61 
to get a numerical approximation of ip m ,n)- in the same graph, we add show the 
coefficient c 2 of the Vicsek model [2UJE2] for d = % ■ The relative error between 
the two curves is very small (around 5%). This similarity between the two curves 
shows a strong connexion between the PTWA model and the Vicsek model. Work 
is in progress to study the link between the two models more deeply. 
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The generalized collisionnal invariant ^ m , n £VW» computed with finite difference 




Figure 2: Left figure: the generalized collisional invariant ifj m , n for A = 1 and a = 1 
computed using m = 30 and n — 61. Right figure: we compute Cip mjn using a 
finite difference method with A9 = .2 and An = .2. We clearly recover the function 
— sin# (see figure [3] for a more detailed comparison). 



The residue \&p„ hn + sin# 




Figure 3: The residual \jCip m>n + sin^loo estimated on the interval (6, k) G [— 7r, 7r] x 
[—5, 5] for different values of (A, a). i/} m>n is computed as in figure [2] (left) and Cip m ,n 
is computed using a finite difference scheme (with AO = An = .2). The residual 
increases with a and decreases with A. 




Figure 4: The coefficient c 2 in the PTWA model (I3.38P computed for A = 1 and 
different values of a (blue) and the coefficient c 2 in the Vicsek model (green). The 
relative error between the two curves is around 5%. 
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5 Conclusion 

In this work, we have introduced a new Individual-Based Model describing the dis- 
placement of individuals which tend to align with theirs neighbors. This model, 
called 'Persistent Turning Walker model with Alignment' (PTWA), is a combina- 
tion of the phenomenological Vicsek alignment model [JT] with the experimentally 
derived PTW model of fish displacement [23]. We have established the macroscopic 
limit of this model within a hydrodynamic scaling where the radius of interaction 
of the agents is tied to the microscopic scale. The derivation uses a new notion of 
'Generalized Collisional Invariant' developed earlier in [20]. The numerical compu- 
tations of the coefficients involved in this macroscopic model have shown that there 
are important similarities between the PTWA model and the Vicsek model at large 
scale. 

The present work proves that the addition of a local alignment rule in the PTW 
model changes drastically the large-scale dynamics as compared to the PTW model 
without alignment interaction. Indeed, while the PTW model without alignment 
is diffusive at large scales, the PTWA model becomes hyperbolic, of hydrodynamic 
type. As a summary, local alignment generates macroscopic convection. 

In future work, the relation between the PTWA and Vicsek dynamics will be 
further explored, both at the microscopic and macroscopic levels. This ensemble 
of models forms a complex hierarchy. Numerical simulations and comparisons over 
a wide range of parameters will be performed to better understand the relations 
between these models. 

Many questions concerning the derivation of macroscopic models remain open 
in this context. One possible route is to explore what the macroscopic limit of the 
PTWA model becomes when an attraction-repulsion rule is added. More generally, 
it may be possible to classify the different types of Individual-Based Models by 
looking at their corresponding macroscopic limits. Another direction is to quan- 
tify how close the macroscopic model is to the corresponding microscopic model. 
In particular, the question of determining what minimal number of individuals is 
required for the macroscopic description to be valid is of crucial importance. All 
these questions call for deeper numerical studies which will permit to understand 
when the microscopic and macroscopic descriptions are similar and when they are 
not. 



5 Conclusion 



27 



Appendix Al: Proof of lemma 13.41 (ii) (functional 
analytic proof) 

Proof. First, we prove the uniqueness of the solution of (13.251) in E. Indeed, we 
have shown in section [3.3.21 that the null space ker(£) of £ consists of the constant 
functions. Therefore, ker(£)r\E = {0}, which shows the uniqueness of the solutions 
of (153511 in E. 

To prove the existence of a solution of (13.251) . we first consider a slightly modified 
version of equation (13.251) : for a given e > 0, we want to solve 

-eiP + CiP = x- (5.1) 

Thanks to this modification, we have the inequality: 

< exjj - Zty, if, >,= e\1>\l + a 2 \d K ^\l > e^. 

Therefore the operator eld — £ is coercive, so we can apply the theorem of J. L. 
Lions in [21] which gives a weak solution ip £ in E of the problem (15. ip . 

To find a solution of £ip = x, we need to extract a convergent subsequence of 
{ipe}e>o when e goes to zero. The limit will satisfy (I3.25p . Since E is an Hilbert 
space, it remains to prove that the family {ip £ } £ >o is bounded in E. For that, we 
proceed by contradiction. If the family {ip £ } £ is not bounded in E as e tends to 0, 
there exists a subsequence e n such that: 

\We n \n > +00 , E n ► 0. 

To simplify the notations, we use the subscript e for e n in the following. Defining 
the functions: 

Ug = w € (5 - 2) 

with N £ = we have that: 

-eU £ + £U £ = ^-. 



N, 



e 



Since the sequence {U £ } £ is bounded (\U £ \^ = 1), we can extract a weakly convergent 
subsequence (denoted by e once again) such that: 

U £ Uq weakly in L 2 . 

In particular, since N £ +00, we have that £Uo = and therefore by uniqueness 
Uq = 0. This means that U £ converges weakly to zero. We will obtain a contradiction 
with the fact \U e \n = 1 if we prove that U £ converges strongly to zero. 
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To prove the strong convergence of U £ , we decompose the functions U £ in two 
parts. For that, we introduce the vector space L: 



It is easy to see that L C E. We denote by L the orthogonal space of L such that: 

E = L®L L . 

We can decompose the sequence U e as U £ = $ e + v £ with $ e G L and t> £ G L" 1 . 
First, we are going to prove that v £ converges to zero using that £ is coercive on 
L L . Taking the scalar product of the equation (15.21) against U £ , we find: 



Since we have the equality < CU £ , U £ > M = — a 2 |<9 K £4| 2 (13.141) and d K U £ = d K v e , we 
obtain that: 



L = {$ e L 2 ^ 1 ) I [ $>{6)M{9) d9 = 0} 



-e\U £ \l+ < CU £ , U £ > M = — < X ,U £ >. 



Therefore, at the limit e — > 0, we have: 



< cu e , u £ > 




0. 




(5.3) 



Then we use the Poincare inequality for Gaussian measures [26J: 




(5.4) 



with C a positive constant and / the mean of / defined as: 




Applying the Poincare inequality ( 15.41) to v £ leads to: 




(5.5) 



Since « £ 6 L 1 , for all $(0) G L, 



we have 
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Therefore v £ {6) = 0. Combining the inequality (15. 5p with (15.31) yields: 



|2 e-> 



, 0. 



It remains to prove that $ e converges to zero. With this aim, we take the scalar 
product of the equation (15. 2p against the function k. Once we take the limit e — » 0, 
we find: 



< CU £ , k 0. 



Using that |<9 K u £ |^ also converges to zero, we deduce that: 

/ K 2 d e U £ fid6dK £ -^0. (5.6) 

We would like to use once again a Poincare inequality. With this aim, we define 
the function h £ {6) as: 

h £ (6) = [ K 2 U £ (6,K)Af(K)dK 



and we use the notation: 

\h(9)\ 2 M = f \h(6)\ 2 Md6. 

J 6 

So equation (15. 6|) can be read as Idgh^ 2 ^ 0. The usual Poincare inequality 
gives: 

\h £ -h £ \ 2 M <C\d e h £ \ 2 M , (5.7) 

with h £ = J g h £ (6)M.(6) d9. But since we already know that U £ converges weakly to 
zero, we have: 

K =< u e , k 2 >^ o. 

Therefore the Poincare inequality (15. 7\i yields h £ — ^ 0, or in other words: 

/ K 2 U £ fid6dK e -^0. (5.8) 
Je,K 

Since v £ converges to zero, equation (15. 8p leads to: 

<$> £ (6)k 2 M(6)M(k) d6dK ^ 0, 

which finally gives that $ £ also converges strongly to zero in L 2 ^. 
Since both v £ and $ e convergence strongly to zero, U £ converges strongly to zero 
as well. This contradicts that \U £ \^ = 1 for all e. Therefore, the sequence ip £ is 
bounded in L 2 ^, so we can extract a subsequence which converges weakly to ip Q in 
L 2 ^. This function ip has to satisfy: 

£ipo = X 

which ends the proof of the lemma. 

□ 
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Appendix A2: Proof of lemma 13.41 (ii) (probabilis- 
tic proof) 

Proof. The operator C is the infinitesimal generator of the following stochastic 
differential equation: 

d6 = K dt, (5.9) 
dK = -A(sin0 + k) dt + \p2adB u (5.10) 

For any function (p regular enough, we can define the semi-group: 

P t (<p)(6,K) = E[<p(X t )\X = (6,k)], 

with X t the stochastic process solution of fl5.9p - fl5.10p . This defines a solution of 
the following equation (see [35]): 

J d t u = Cu 

\ U t=0 = if. 

In particular, if we define u(t) = Pt(x)> a simple integration by part leads to: 

u(t)- X = I Cu(s)ds. (5.11) 
Jo 

Therefore, we will find a solution to (13.251) if we are able to prove that u{t) 0. 
For that, we first notice that the equilibrium measure associated with £ is given by 
ji (I3.13P and its adjoint operator in is given by:. 

£*ip = —ndeil' + A sin dd K i^ — \Kd K ijj + a 2 d^tjj. 

Moreover, we can find a Lyapunov function associated with C. The function 
V(9, k) = 1 + k 2 satisfies: 

C*V = 2\sm9K-2\K, 2 + 2a 2 

< 2Xk - A(l + k 2 ) - Xk 2 + A + 2a 2 

< -XV + 2(a 2 + X)t { ^ 2+VIT ^^ } . 

Therefore V is a Lyapunov function in the sense of [TJ Def. 1.1]. Since B = 



S 1 x{|ft;|<2 + yl + (2a /A + l) 2 } is compact, B is a "petite set" in the terminology 
[H Def. 1.1] of Meyn & Tweedie [32]- So we can apply [H Th. 2.1] and conclude 
that there exists a constant K2 > such that for all bounded function tp satisfying 
J g K <p // d6dn = 0, we have: 



\Pt(<p)\* < K 2 \\ 



-XI 
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Therefore, we can pass to the limit t — > oo in ( 15.111) to find that: 



-X = I £u(s) ds, 
Jo 

Defining the function if> — — / °° u(s) ds, we get a solution to: 

= X- 

For the uniqueness of the solution, we proceed as in appendix Al. 



□ 



Appendix B: Proof of proposition SUL 

Proof. We first prove the following lemma. 

Lemma 5.1 For every integer m and every positive integer n > 0, we have: 



C(cp m P n )= J2 D m > n (j,k)<p m+j Pn+k 

with D m ' n a 3 x 3 matrix given by: 



-i<i<i 

-Kfc<l 



D' 








fix my/n 


—An 








h Vn + 1 



with: 



R - %a 



-f3 2 Vn+1 
4a 



(5.12) 



Proof. First, using the properties of the Hermite polynomial^!, we can find several 
properties of P n : 



K = —y/nP n -i, 
a 

nP n = -j= [y/n + lP n +i + VnP n - 



(5.13) 



2 Indeed H' = nH n -i and xH n = H n+ i + nH n -i 
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In particular, the polynomials P n are eigenfunctions of the self-adjoint part of C: 

- \nd K P n + c?d\P n = -XnP n . (5.14) 

Then, we compute: 

C{ip m P n ) = KP n d e ip m A sin 9 ip m d K P n + ip m ( Xk3 k P n + a 2 d 2 k P n ). 
The derivative of ip m with respect to 9 is given by: 

irnO \ / „irn6 \ „im9 I i A 



^2 gim0 /g*^ e _ *^ 1 \ 



1 -^sin0.M 



2a 2 V 2i V^M 



= imif m - 4^(</Wi - <fm-l)- 
Using (15.131) . we also have: 

KP n d e (f m = (y/n + lP n +i + y/nP n -i) (im<p m - ^(vWi - <^ m -i)j • 

= 7^ (mVn + lPn+i^m + mv/nPn-i^m) (5.15) 

(Vn + l-Pn+l^m+l + VnPn-l^Pm+l) (5.16) 

(A/n+IP n+ iV9 m _i + v/riPn-i^m-i) • (5-17) 



Thus, we have: 



A sin 9 (p m d K P n = -A — ip m VnP n -i 

1% a 

= \ /n{(p m +l - (p m -l)Pn-l- (5.18) 

la 

Finally, since P n satisfies (I5.14p . we get: 

ip m (-\nd K P n + a 2 d 2 k P n ) = -\n ip m P n . (5.19) 

Combining flBTTol) (15TTBD (IBTT7D (1BTT8D (15TT91 . we find the expression (I5^2]l of D m ' n . 
□ 
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To find the matrix representation of the operator C in V m>n , we introduce the vectors 
u and v defined by: 



u 

v 



(Po, 



) Pn 



With these notations, a function ip e V mjn with a matrix representation X (14. 4p 
can be written as: 

^m,n = X! Cj(pjP k = U T XV. 

\j\<m,0<k<n 

Moreover, thanks to the matrices defined in (14.5)) . we can write for example 

' V ' • E C'ifjl'k I 

\j\<m, l<k<n 

E jcfcp*. 



u 1 XL-iV 
(D lU ) T Xv 



\j\<m, l<fc<n 

For a function ip G L^, using the lemma I5TT1 we can write: 



EC 



/ 



\ 



-1<?'<1 
-Kfe<l 



= EC (a 

m,n 

+/3 2 (y?m-i - ¥Wi)\/ra + lP n+ i - \ip m nP n ] 
Therefore, for every ip e V^j i7l , we have: 

<£^> M = < ^(D.ufXi^L^+P^ufXiL^y/D^v) 

-/3 2 (L_ 1 u) T X( A /^L_ 1 u) + p 2 {L +1 u) T X{jD~ 2 L_ lV ) 

+f3 2 (L_ lU fX(L +1 ^D 2 v) - f3 2 (L +1 u) T X(L +1 JF 2 v) 
-\u T X(D 2 v) , > At . 

We can simplify this expression: 

<C4>, ip> fl = < u T p 1 D 1 X{^ 2 L- 1 + L +1 ^ 2 )v 
+u T /3 2 (-L_i + L + i) T X(/5^L_i) v 
+u T /9 2 (L_! - L +1 ) T X(L +1 JI)~ 2 ) v 



-\u T X(D 2 v) 
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which finally gives 

< A/>, V >„=< PiMxXNi + f3 2 M 2 XN 2 - \XD 2 , <p >„, 

with Mi, M 2 , Ni and defined in ( 14. 7p . Therefore, using tp = — sin 6, we find 
that X has to satisfy equation (14. 6p . 
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